Synchronization of extended chaotic systems with long-range interactions: 
an analogy to Levy-flight spreading of epidemics 
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Spatially extended chaotic systems with power- law decaying interactions are considered. Two 
coupled replicas of such systems synchronize to a common spatio-temporal chaotic state above a 
certain coupling strength. The synchronization transition is studied as a nonequilibrium phase 
transition and its critical properties are analyzed at varying the interaction range. The transition 
is found to be always continuous, while the critical indexes vary with continuity with the power 
law exponent characterizing the interaction. Strong numerical evidences indicate that the transition 
belongs to the anomalous directed percolation family of universality classes found for Levy-flight 
spreading of epidemic processes. 
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Since its discovery Q , chaotic synchronization has be- 
come a very active and important field of research Q, 
especially for its applications in control and secure com- 
munications 0, B • It has been observed in various ex- 
perimental setups from semiconductor lasers to chemi- 
cal reactions 5j. The synchronization transition (ST) is 
particularly interesting in spatially extended chaotic sys- 
tems 6] where, due to the spatial dependence, a parallel 
with nonequilibrium phase transitions can be drawn. Ex- 
tensive numerical studies of coupled map lattices (CML), 
which are prototype models for space-time chaos 7] , have 
shown that ST between two replicas of CML with short 
range interactions is continuous and the associated criti- 
cal properties belong to two distinct universality classes 
depending on the nature of the individual elements of 
the system H, H, [lOj El • For smooth chaotic maps, lin- 
ear mechanisms rule information propagation, and ST oc- 
curs as soon as the transverse Lyapunov exponent crosses 
zero. In this case, the transition belongs to the mul- 
tiplicative noise universality class 12]. For discontinu- 
ous maps, a linear stability analysis is no more sufficient 
to locate the ST, and the replicas synchronize (for defi- 
nitely negative transverse Lyapunov exponent) when the 
spreading velocity of finite amplitude perturbations van- 
ishes gj. The synchronization in such systems is con- 
trolled by nonlinear effects and the transition belongs 
to the directed percolation (DP) universality class Pjj, 
analogously to what found for cellular automata Eil . At 
present we thus have a satisfactory understanding of the 
critical properties of ST for short-range coupled chaotic 
systems |9l llflj. 

The aim of this Letter is to study the synchronization 
transition of extended chaotic systems with long range 
spatial coupling. Long range interactions are indeed rel- 
evant to many real systems such as neuron populations 



E^, Josephson junctions 16], cardiac pacemaker cells 
[13 , and to the issue of disease spread triggered by avia- 
tion traffic We consider, as a prototype of systems 
with long range interacting elements, a model of CMLs 
with coupling decaying as a power-law 

ei Him. As 

a first step we investigate the case of maps whose dy- 
namics is dominated by nonlinear mechanisms. We show 
that the critical indexes characterizing the transition de- 
pend continuously upon the exponent a, which sets the 
coupling range. Moreover, strong numerical evidences 
indicate that the critical exponents are compatible with 
those found in a model for anomalous directed percola- 
tion studied by Hinrichsen and Howard [22j (HH in the 
following) . 

Anomalous DP arises in epidemic spreading when- 
ever the infective agent can perform unrestricted Levy- 
flight sJS EH • Such processes, originally introduced in 
Ref. 24], can be modeled by assuming, e.g. in d— 1, 
that an already infected site propagates the disease to 
any other site with a probability P(r) ~ r -{i+v) alge- 
braically decaying with the distance r (similarly to the 
distance distribution of human travels, which are nowa- 
days responsible for the geographical spread of infectious 
diseases [HI). Numerical studies of a stochastic lattice 
model, which generalizes directed bond percolation, have 
shown that the critical exponents vary continuously with 
g |22|, despite the divergence of the average interaction 
distance for a < 1 and of the second moment of P(r) 
for a < 2. These findings confirm previous theoreti- 
cal results [23| according to which usual DP should be 
recovered for sufficiently short-range coupling (namely 
for a > a c = 2.0677(2)) and the mean-field description 
should become exact for a < a m = 0.5. Furthermore, 
close to a m , the critical exponents found by HH agree 
with renormalization group calculations [2^ . 
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We investigate the synchronization transition of two 
coupled replicas of CML defined as: 

Xi(t+1) = (l-7)F(xi(t))+<yF(yi(t)) 
2/*(* + 1) = {l-l)F{Ut))+lF(xi{t)) (1) 
where £, = {x~i, y,} is given by: 

M 

The labels t and i are the discrete temporal and spatial 
indexes, Xi(t),yi(t) G [0 : 1] are the state variables, L is 
the lattice size (i = 1, . . . , L), and periodic boundary con- 
ditions are assumed (i.e. Zi + L{t) = Zi(t)). As local dy- 
namics, throughout this work, we employ the Bernoulli 
map F(x) = 2x (mod 1), which is characterized by a 
discontinuity that makes the nonlinear effects dominat- 
ing |23|. The interaction range is controlled by the ex- 
ponent a: for <j = globally coupled maps are obtained, 
while standard CML's with nearest neighbor coupling are 
recovered in the limit cr — » oo. Two coupling constants 
e and 7 enter the system definition. The former mea- 
sures the amplitude of the spatial coupling and is fixed 
to 2/3. The latter controls the strength of the site- wise 
interaction between the two replicas. In the range of cr 
examined, a synchronized state is always achieved as the 
coupling exceeds a certain >y c = >y c (q, cr). 

Usuall y th e model has been studied in the fully coupled 
case EMEU, i-e. for j m (q) = m with M = (L-l)/2. 
Instead by following Ref. 26], we consider here a modi- 
fied version of the coupling, where j m (q) = q m — 1 with 
q- values typically chosen as q = 2,4 and 8. The sum 
in Eq. (0 extends to M = log (L/2); finally rj(a) = 
^Em=iM(imfe)) _(J ^ S a normalization factor. The ra- 
tionale for the choice j m (q) =q m —l lies in its convenience 
from a computational point of view. As a matter of fact, 
the chosen coupling scheme allows for simulating very 
large systems: each updating step can be performed in 
0(L \og q L) operations instead of 0(L 2 ) needed in the 
fully coupled case. Moreover, a simple reasoning shows 
that the model with j m (q) =q m — l and exponent a ex- 
hibits the same critical properties of the fully coupled one 
with a different exponent crf c = a + 1. Indeed both ver- 
sions of the model should display the same critical prop- 
erties once the spatial interactions scale analogously. For 
the modified model the coupling weight over the interval 
]jm(o)^ jm+i(o)] containing a single coupled site is simply 
l/j m+ i(g) cr ~ q-{m+i)(j ^ w hii e f or f u r}y coupled case 
this amounts to k ~ afc ~ q -(™+iK* fc -i) , 

Therefore, the two weights scale in the same manner only 
for <t/ c = <t+1 as confirmed by numerical results reported 
in the following. 

To characterize the ST it is useful to introduce the 
synchronization error Wi(t) = \xi(t) — y%(t)\. Its spatio- 
temporal evolution is reported in Fig.Qfor various values 




FIG. 1: (Color online) Space-time plot of synchronization 
error Wi(t) = \xi(t) — y%(t)\ for two replicas of the CML (1-2) 
for q — 2 and M — 9 with 7 slightly larger than j c and differ- 
ent cr's, see labels. The space (resp. time) axis correspond to 
the horizontal (resp. vertical) axis. The colors code the ab- 
solute value of the synchronization error in logarithmic scale. 
Colors go from intense white (yellow online) corresponding to 
Wi(t) rsj 1 to dark black associated to the fully synchronized 
state, i.e. Wi(t) = 0. 

of a just above the synchronization transition (7 > j c ). 
In the short range limit a — > 00, percolating structures, 
typical of ordinary DP, are clearly observable. As the 
range of the interaction increases (i.e. a decreases) the 
spatial structures tend to be smoothed out. Finally, for 
cr < 0.5, the spatial structures are no more detectable. 

To be more quantitative, we consider the spatial av- 
erage of the synchronization error p 1 (t) = Wi (t) / L 
which is the natural order parameter: it vanishes, at suf- 
ficiently long times, whenever a complete synchronization 
is achieved, i.e., p* = lim^oo p 1 {t) = for 7 > 7 C ; whilst 
it remains finite at any time in the unsynchronized state, 
i.e. for 7 < 7 C . In terms of p 7 (t), we can now define the 
critical exponents 5, (5 and z that will be function of a: 5 
rules the temporal scaling of the order parameter at the 
critical point 7 = j c i.e. p lc (t) °° P controls the way 
the asymptotic value of the order parameter approaches 
zero, p* ~ (7 C — 7)^ for 7 < 7 C ; and finally the dynami- 
cal exponent z can be defined in terms of the finite-size 
scaling relation, valid at the critical point: 

p ye (t)~L-*'f(t/L'). (3) 

These exponents, which fully characterize the transi- 
tion (ruling the divergence of, e.g., the spatial and tem- 
poral correlation lengths at the critical point), have been 
measured for different values of cr, by means of the follow- 
ing procedure (see Fig. [2t,b and c for an illustration of 
the method and to appreciate the quality of the results). 
The critical coupling 7 C is estimated by examining the 
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FIG. 2: (Color online) (a) Average synchronization error p 1 
vs time for the CML with q = 4 and M = 8 for a = 1.4, the 
curves refer to 7 = 0.2977, 0.29789, 0.29792, 0.29797, 0.2981 
(from top to bottom). Here j c = 0.29792(2). The straight 
line indicate the slope S = 0.33(2) as obtained by a best fit. 
(b) Asymptotic synchronization error p* vs j c — 7 for the 
same parameters as in (a), the straight line indicates the slope 
p = 0.49(5). (c) Finite-size scaling p lc (t)L Sz vs t/L z from 
which we estimate 2 = 1.05(5), symbols refer to different sizes 
L = 2 q M (here q = 4) and M as in the legend, (d) The 
synchronization error p lc vs time at criticality for a — 2 using 
q = 2 with M — 16 and q = 4 with M — 8 together with the 
results for the fully coupled model for af c = 1 + cr = 3 and size 
L = 8, 193. The top straight-line displays the slope 5 = 0.2 
as estimated by a best fit; the bottom straight-line shows 
the (short range) DP value 5 = 0.15945. Data are obtained by 
averaging p 1 {t) over many (~ 100) different initial conditions. 



temporal behavior of p 7 (£), averaged over many random 
initial conditions, for different values of 7. In particular, 
7 C is determined (typically with an absolute error within 
the range [2 x 10 -5 : 2 x 10 -4 ] ) as the point at which 
a scaling behavior for p 1 (t) emerges (see Fig. Et)- With 
this method 8 is directly computed by means of a best 
fit procedure on p lc (t). Once 7 C is known, (3 is obtained 
by measuring the saturation values p* in the desynchro- 
nized regime (i.e. for 7 < j c ) for sufficiently long chains, 
see Fig. Efc. Finally, performing several simulations with 
different sizes at the critical point 7 = 7 C , the exponent 
z is estimated through a data collapse based on Eq. JSJ. 
Although this procedure works satisfactorily (Fig. it 
is affected by large errors due to the unavoidable subjec- 
tivity in judging the collapse. 

We remark that, as typical of strongly nonlinear 
maps H 0, 7 C cannot be identified by the condition 
of vanishing transverse Lyapunov exponent and, more- 
over, it depends not only on a but also on q. However, 
the critical exponents, apart from finite-size effects, are 
not influenced by q as one can judge from Fig. 01, where 
the synchronization error at the critical point is shown 
for q = 2,4 and a = 2. In the figure we also display 



FIG. 3: (Color online) Critical exponent /3 vs a obtained from 
best fits using different lengths and basis q = 2, 4, 8. The black 
thick curve represents the data of HH 22] , the bottom straight 
line is the DP value, while the top one represents the mean 
field prediction by HH. CML with sizes from L = 2 16 « 6x 10 4 
up to L = 2 x 8 7 « 4 x 10 6 maps have been used. Inset: 
Critical exponent 8 vs a symbols are as in the main panel. 

the result for the fully coupled model with the exponent 
ruling the interaction range set to <jf c = 3. The data 
confirm that cr/ c = 1 + cr as previously discussed. 

In Fig. 03 we display the exponents (3 and 8 as a func- 
tion of cr, and compare them with the results of HH 
obtained for stochastic models of Levy-flight epidemic 
spreading [2^. The agreement between the two sets of 
data is rather good. Moreover, by increasing q and so 
the maximal affordable spatial length L, we observe that 
the exponents approach those of HH. This leads us to 
conclude that discrepancies are mainly due to finite-size 
effects. The results for z (not reported) also confirm the 
overall agreement. Due to the strong finite-size effects 
and to computational limitations it was not possible to 
analyze systems with a < 0.5. For anomalous DP, it has 
been shown that in the interaction range [<j m : a c ] the 
following scaling relation holds EH 

A = l-cr+(l-2£)z = 0. (4) 

As shown in Table [I] the relation (0J is fulfilled within the 
errors, strengthening the parallel between ST in spatio- 
temporal chaotic systems with power-law interactions 
and anomalous DP. It is worth stressing that it is abso- 
lutely nontrivial that the deterministic system here inves- 
tigated was in quantitative agreement with the stochastic 
one studied in HH. The fact that the critical exponents 



° [q] 


3.0 [2] 


2.0 [4] 


1.4 [4] 


1.2 [4] 


0.8 [8] 


0.5 [8] 


A 


0.93(6) 


-0.13(9) 


-0.04(12) 


-0.03(10) 


-0.03(23) 


0.04(13) 



TABLE I: Scaling relation Q) for various values of a, for 
each measurement the corresponding basis q is reported. The 
measured A are compatible with zero in the range I cr m ■ I • 
See text for details. 
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are the same suggests that the universality class (actu- 
ally the family of universality classes labeled by a) is the 
same for both models. 

The agreement of the critical exponents with those of 
Ref. 22] (shown in Fig.© represents a strong indication 
that the synchronization dynamics associated to model 
(|H2j) should admit the same stochastic field description 
proposed for anomalous directed percolation [2^, |2J • 

d t n = D N V 2 n + L> A V a n + rn - An 2 + g(ri)£ (5) 

where, in the present context, n represents the coarse 
grained synchronization error obtained by averaging 
Wi(t) over a suitable space-time cell. The parameters en- 
tering Eq. © are the standard and anomalous diffusion 
coefficients Dn and Da, respectively; r that measures 
the distance from the critical point (i.e. from the syn- 
chronization threshold j c ) and the amplitude of the non- 
linear term A. The fractional derivative V a is non-local 
and can be defined through its action in Fourier space 
— \k\ a . Finally, £ is a zero- average ^-correlated (in 
space and time) Gaussian noise field with unit variance 
and g 2 (n) oc n. For systems with short range interac- 
tions, Eq. © with Da = (i.e. the field description 
of ordinary DP) was heuristically derived in from a 
simple stochastic model, mimicking the behavior of Wi(t) 
for Bernoulli maps. That derivation can be straightfor- 
wardly extended to the present case giving Eq. ® as 
a macroscopic description of the dynamics. Indeed the 
only difference comes from the spatial coupling that can 
be handled as in (2^ (where the synchronization of long 
range coupled chaotic oscillators was considered) leading 
to the fractional derivative term. However, a detailed 
derivation goes beyond the scopes of the present Letter. 

In conclusion, we have provided strong numerical ev- 
idence that the synchronization transition of two cou- 
pled replicas of a chain of discontinuous maps with long 
range interactions is characterized by a continuum of uni- 
versality classes labeled by the exponent a, which sets 
the interaction range. Moreover, the critical exponents 
agree with those reported for anomalous directed per- 
colation [2^. Preliminary results indicate that also for 
continuous maps the critical exponents depends on a and 
that they do not belong to the anomalous DP universal- 
ity class. By analogy with the present results, we con- 
jecture that the critical properties of the synchronization 
transition for continuous maps should be reproducible in 
terms of a multiplicative noise Langevin equation with 
long-range interactions (i.e. an equation similar to (© 
but with g(n) oc n [13)- To the best of our knowledge 
this equation has not been yet subject of theoretical or 
numerical analysis. It will be a definitely worth task to 
undertake such study in the future 
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